\appendix
\chapter{Wizualizacja kształtu powierzchni}
\begin{quote}
\begin{flushright}
Żadnym pomiarom nie można wierzyć,\\ dopóki nie zostaną potwierdzone przez teorię.\\
SIR ARTHUR EDDINGTON
\end{flushright}
\end{quote}


Sprawą o zasadniczym znaczeniu podczas modelowania jest możliwość dokonania oceny kształtu otrzymanej krzywej lub powierzchni. Ocena ta jest podstawą do zaakceptowania projektu lub do wprowadzenia poprawek. Najprostszym środkiem do dokonania takiej oceny powierzchni jej "realistyczny" obraz, wykonany z użyciem dowolnego modelu oświetlenia. Model taki symuluje określone własności optyczne powierzchni i stosownie do potrzeb umożliwia wybieranie źródeł światła. Jednak w praktyce wygląd powierzchni na otrzymanych w ten sposób obrazach często nie wystarczy do jej oceny.

Ocena kształtu powierzchni może obejmować trzy podstawowe elementy. Pierwszy z nich to rozkład krzywizny na powierzchni, a w szczególności obecność zafalowań. Drugi to obecność nieciągłości płaszczyzny stycznej lub krzywizny powierzchni na połączeniach płatów wchodzących w jej skład. Przedstawione niżej metody wizualizacji kształtu mogą być stosowane do oceny obu tych elementów, choć w każdym z tych przypadków może być potrzebne utworzenie nieco innych obrazów. Trzecim elementem, czyli utrzymaniem tolerancji odtworzenia teoretycznego kształtu powierzchni elementu mechanizmu współpracującego z innym elementem (na przykład w taki sposób, że powierzchnie tych elementów ślizgają się po sobie), nie będziemy się tu zajmowali.

Istnieją dwa podstawowe sposoby obrazowania kształtu powierzchni. Pierwszy to określenie pewnej \textbf{funkcji kształtu} $f$ na powierzchni i przedstawienie tej funkcji za pomocą kolorów. Jeśli odwzorowanie wartości funkcji $f$ na kolory, zwane \textbf{paletą}, jest ciągłe, to otrzymany obraz uwidocznia nieciągłości funkcji $f$. Nieciągła paleta umożliwia otrzymanie obrazu warstwie funkcji $f$. Drugi sposób wizualizacji polega na narysowaniu obrazu \textbf{krzywych charakterystycznych} leżących na powierzchni. Krzywe te w szczególności mogą być określone za pomocą pewnej funkcji kształtu, na przykład mogą stanowić jej warstwice albo linie najszybszego spadku.

\textbf{Obrazem funkcji kształtu} nazwiemy rysunek powierzchni (w pewnym rzucie) pokolorowany w ten sposób, że obraz punktu $\boldsymbol{p}(u,v)$ ma kolor z palety odpowiadający wartości $f(u,v)$. Aby uniezależnić obraz, na podstawie którego będzie dokonywana ocena kształtu, od parametryzacji płata, wystarczy określić funkcję $f$ w ten sposób, że dla każdego punktu $(u,v) \in A$ wartość $f(u,v)$ jest określona przez kształt płata w otoczeniu punktu $\boldsymbol{p}(u,v)$. Istnieją różne funkcje spełniające ten warunek, na przykład krzywizny powierzchni. 

Właściwie dobrane narzędzie wizualizacji, czyli odpowiednia funkcja kształtu lub rodzina krzywych charakterystycznych, uwypukla na obrazie negatywne cechy kształtu powierzchni. Interpretacja takich obrazów bywa trudna i wymaga od projektanta dużego doświadczenia i wprawy. Prawdopodobnie największą sztuką jest podjęcie decyzji, jak zmodyfikować powierzchnię, aby skutek tej modyfikacji był zadowalający.

\section{Funkcje kształtu i ich warstwice}
\subsection{Własności warstwic}
Rozważmy płat $\boldsymbol{p}: A \to \mathbb{E}^3$ i krzywą $\boldsymbol{b}:(a,b) \to A$ w dziedzinie płata $\boldsymbol{p}$. Krzywa $\boldsymbol{c}(t) = \boldsymbol{p}(\boldsymbol{b}(t))$ leży na płacie $\boldsymbol{p}$ (zobacz p. C.2.3). Przypuśćmy, że parametryzacja płata $\boldsymbol{p}$ jest klasy ${C}^n$. Jeśli krzywa $\boldsymbol{b}$ jest klasy ${C}^k$, to oczywiście rząd ciągłości jest nie mniejszy niż $min(n,k)$ (uwaga: nie wyklucza to ciągłości wyższego rzędu krzywej $\boldsymbol{c}$). O klasie ciągłości parametryzacji płata $\boldsymbol{p}$ można by więc wnioskować na podstawie klasy ciągłości krzywych leżących na tej powierzchni. Nie jest to praktyczna metoda, ponieważ obraz krzywej nie daje informacji o jej parametryzacji, a ponadto gdy celem jest badanie kształtu, należy je uniezależnić od parametryzacji używanej w danej chwili do reprezentowania płata. Zasada badania jest jednak podobna po sformułowaniu problemu za pomocą pojęcia ciągłości geometrycznej (zobacz rozdziały 8 i 9), które jest niezależne od parametryzacji.

Załóżmy, że płat $\boldsymbol{p}$ jest klasy $G^n$ dla pewnego $n \geq 1$, a zatem istnieje jego regularna parametryzacja klasy $C^n$. Jeśli nie istnieje regularna parametryzacja klasy $C^{n+1}$ (co oznacza, że powierzchnia nie jest klasy $G^{n+1}$), to dowolną regularną parametryzację klasy $C^n$ nazwiemy \textbf{parametryzacją maksymalnie gładką}. W dalszym ciągu będziemy rozważać funkcje kształtu określone w dziedzinie powierzchni o takiej właśnie parametryzacji. Niech $f: A \to \mathbb{R}$ oznacza funkcję klasy $C^k$ ($k \geq 1$), która ma niepusty zbiór miejsc zerowych w zbiorze $A$ i dla każdego $(u,v) \in A$, takiego że $f(u,v) = 0$, jest spełniony warunek $f_{u}(u,v) \neq 0$ lub $f_{v}(u,v) \neq 0$. Z \textbf{twierdzenia o funkcji uwikłanej} wynika, że jeśli $f(u_{0},v_{0}) = 0$, to istnieje przedział $(a,b)$ i regularna krzywa $\boldsymbol{b}:(a,b) \to A$ klasy $C^k$, która przechodzi przez punkt $(u_{0},v_{0})$ i której wszystkie punkty są miejscami zerowymi funkcji $f$. Ponieważ zarówno krzywa $\boldsymbol{b}$, jak i płat $\boldsymbol{p}$ są regularne, więc krzywa $\boldsymbol{c}:\boldsymbol{c}(t) = \boldsymbol{p}(\boldsymbol{b}(t))$ klasy $C^m$ ($m$ = min$(n,k)$) jest regularna, a zatem jest klasy $G^m$.

Przypuśćmy, że funkcja kształtu jest określona w taki sposób, że jeśli powierzchnia jest klasy $G^n$ (dla dowolnego $n$), to funkcja kształtu $f$ jest klasy $C^{n-l}$ (dla pewnego $l$ niezależnego od $n$). Łatwo możemy sprawdzić, że jeśli krzywa $\boldsymbol{c}$, określona przez zbiór miejsc zerowych funkcji $f$, jest klasy $G^m$, ale nie $G^{m+1}$, to powierzchnia $\boldsymbol{p}$ jest klasy co najwyżej $G^{m+1}$ (uwaga: należy jednak upewnić się, że w żadnym punkcie krzywej funkcja $f$ nie ma obu pochodnych cząstkowych równych 0).

Mając ustaloną funkcję $f$, możemy określić w dziedzinie płata rodzinę krzywych wzorem $f(u,v) - c = 0$, w którym $c$ jest parametrem. Krzywe te są \textbf{warstwicami} funkcji $f$. Krzywe powstałe w wyniku odwzorowania warstwic na powierzchnię są często używane do badania kształtu. Ich własności możemy oczywiście badać, korzystając z twierdzenia o funkcji uwikłanej, tak samo jak obraz szczególnej warstwicy odpowiadającej wartości $c = 0$. W szczególności wnioski otrzymane wyżej dla krzywej $f(u,v) = 0$ pozostają prawdziwe dla dowolnej warstwicy funkcji $f$.

Gładkość krzywych charakterystycznych, czyli klasa ich ciągłości geometrycznej, jest tylko jednym z elementów oceny. Inne elementy to kształt (np. zafalowania) tych krzywych i rozkład krzywych na powierzchni, a także różnice kształtu poszczególnych krzywych należących do danej rodziny. Sposób otrzymywania obrazów warstwic funkcji kształtu jest opisany w p. F.2.1 i F.2.2.

\subsection{Przekroje powierzchni}
Niech $\boldsymbol{v \neq 0}$ oznacza ustalony wektor, a $\boldsymbol{o}$ punkt (na przykład początek układu współrzędnych). Funkcja kształtu określona wzorem 
\begin{equation}
f(u,v) = \langle \boldsymbol{v}, \boldsymbol{p}(u,v) - \boldsymbol{o} \rangle,
\end{equation}
ma obraz niezależny od parametryzacji płata. Jest tak dlatego, gdyż wartość funkcji $f$ dla ustalonego punktu płata $\boldsymbol{p}(u,v)$ określa jedną ze współrzędnych tego punktu (w pewnym układzie współrzędnych, którego jedna z osi ma kierunek wektora $\boldsymbol{v}$). 

Wektor $\boldsymbol{v}$ jest wektorem normalnym pewnej rodziny płaszczyzn równoległych. Dowolna warstwica funkcji $f$ odpowiada krzywej leżącej w jednej z płaszczyzn tej rodziny. Projektant może dowolnie wybierać kierunek wektora $\boldsymbol{v}$, tak aby otrzymać przekroje, które najlepiej uwidoczniają interesujące go szczegóły kształtu powierzchni. 

Jeśli powierzchnia jest klasy $G^n$, to w każdym jej punkcie, w którym wektor $\boldsymbol{v}$ \textit{nie jest} jej wektorem normalnym, przekrój płaski tej powierzchni jest krzywą klasy $G^n$.

[rysunek]

\subsection{Lambertowskie odbicie światła i izofoty}
Rozważmy najprostszy model odbicia światła od powierzchni, badany przez J. H. Lamberta w 1760 r. W modelu tym jasność $I$ punktu na powierzchni wyraża się wzorem\\
\begin{equation*}
I = c(I_{a} + I_{d}max(0,\cos\measuredangle(\boldsymbol{n},\boldsymbol{l}))),
\end{equation*} 
w którym $I_{a}$ oznacza tzw. intensywność światła rozproszonego w otoczeniu (o której dla uproszczenia często zakłada się, że jest taka sama w każdym miejscu), wektor $\boldsymbol{l}$ określa kierunek, z którego dochodzi światło o intensywności $I_{d}$, wektor $\boldsymbol{n}$ jest wektorem normalnym powierzchni, a współczynnik $c \in [0,1]$ określa, jaka część energii światła padającego na powierzchnię ulega odbiciu. Lambertowskie odbicie światła jest zatem idealnie rozproszone i jasność punktu na oświetlonej powierzchni nie zależy od położenia obserwatora. W szczególności nie ma w tym modelu odblasków spowodowanych zwierciadlanym odbiciem światła.

Jasność punktów powierzchni nie zależy od parametryzacji, a zatem spełnia warunki, jakie powinna spełniać funkcja kształtu. Jej warstwice są zbiorami punktów, w których kąt między wektorem normalnym a ustalonym wektorem $\boldsymbol{l \neq 0}$ jest stały. Jeśli kierunek padania światła na wszystkie punkty powierzchni jest taki sam (źródło światła znajduje się w dużej odległości), to warstwice intensywności światła odbitego od powierzchni, tak zwane \textbf{izofoty}, w rozważanym modelu odbicia światła są warstwicami funkcji kształtu
\begin{equation}
f(u,v) = \langle \frac{\boldsymbol{n}(u,v)}{||\boldsymbol{n}(u,v)||_{2}},\boldsymbol{l} \rangle
\end{equation}
Warstwice tej funkcji dają użyteczną informację o kształcie powierzchni także na jej nieoświetlonych częściach. Funkcja $f$, określona dla powierzchni klasy $G^n$ (o regularnej parametryzacji klasy $C^n$), jest klasy $C^{n-1}$, co wynika stąd, że funkcja $\mathbf{n}(u,v) = \boldsymbol{p}_{u}(u,v) \wedge  \boldsymbol{p}_{v}(u,v)$, opisująca wektor normalny, jest klasy $C^{n-1}$. Zatem izofoty w otoczeniu punktów, w których co najmniej jedna z pochodnych cząstkowych funkcji $f$ jest różna od zera, są krzywymi klasy $G^{n-1}$.\\
(Rysunek F.2)

Przykład funkcji opisującej intensywność oświetlenia powierzchni w modelu odbicia lambertowskiego jest przedstawiony na rysunku F.2. Widzimy na nim powierzchnię klasy $G^1$, pokazaną na rysunkach 9.39 i F.1. Z wyjątkiem miejsc, w których tekstura uwidocznia brzegi płatów wielomianowych, z których składa się powierzchnia, paleta użyta do wykonania obrazka z prawej strony jest ciągła (z dokładnością do technicznych możliwości reprodukcji odcieni). Dzięki zastosowaniu palety użytej z prawej strony kształt izofot, będących brzegami jasnych i ciemnych obszarów, jest wyraźnie widoczny. Na połączeniach płatów, na których krzywizna jest nieciągła, są widoczne nieciągłości kierunku stycznej do izofot.

Przypadkiem szczególnym izofot są tak zwane \textbf{krzywe sylwetkowe}, czyli zbiory punktów, w których wektor normalny $\boldsymbol{n}$ jest prostopadły do ustalonego wektora $\boldsymbol{l}$. Dla powierzchni obserwowanej z dużej odległości (z kierunku $\boldsymbol{l}$) na ogół są one brzegiem widocznej i niewidocznej części powierzchni.
\subsection{Linie odblasku}
Niech $l$ oznacza ustaloną prostą w przestrzeni, określoną za pomocą dowolnego punktu $\boldsymbol{a}$ i wektora $\boldsymbol{t}$. Za jej pomocą określa się krzywą, zwaną linią odblasku. Istnieją dwa sposoby określenia takiej krzywej: zależny od obserwatora i niezależny.

\textbf{Linia odblasku zależna od obserwatora} jest zbiorem punktów powierzchni $\boldsymbol{p}$, w których obserwator (znajdujący się w punkcie $\boldsymbol{e}$) widzi obraz pewnego punktu prostej $l$ w odbiciu zwierciadlanym w tej powierzchni (rys. F.3). Dla dowolnego punktu $\boldsymbol{p}$ na powierzchni określimy następująco wektory:
\\***rys***\\
\begin{equation*}
\boldsymbol{v}=\boldsymbol{e}-\boldsymbol{p},\qquad
\boldsymbol{r}=2
\frac{\langle \boldsymbol{n}, \boldsymbol{v} \rangle}
{\langle \boldsymbol{n}, \boldsymbol{n} \rangle}
\boldsymbol{n}-\boldsymbol{v},\qquad
\boldsymbol{s}=\boldsymbol{a}-\boldsymbol{p}.
\end{equation*}
Wektor $\boldsymbol{v}$ określa kierunek do obserwatora, a wektor $\boldsymbol{r}$ jest jego obrazem w odbiciu symetrycznym względem prostej o kierunku wektora normalnego powierzchni $\boldsymbol{n}$ (na rysunku są zaznaczone tylko kierunki wektorów $\boldsymbol{v}$ i $\boldsymbol{r}$; ich długości są tu nieistotne). Punkt $\boldsymbol{p}$ leży na linii odblasku wtedy, gdy wektory $\boldsymbol{r}$, $\boldsymbol{s}$ i $\boldsymbol{t}$ leżą w jednej płaszczyźnie (czyli są liniowo zależne), a zatem linia odblasku jest zbiorem miejsc zerowych funkcji
\begin{equation}
f(u,v) = det[\boldsymbol{r},\boldsymbol{s},\boldsymbol{t}].
\end{equation}
\textbf{Uwaga:} Liniowa zależność wektorów $\boldsymbol{r}$, $\boldsymbol{s}$ i $\boldsymbol{t}$ jest równoważna liniowej zależności wektorów $\boldsymbol{-r}$, $\boldsymbol{s}$ i $\boldsymbol{t}$. Wektor $\boldsymbol{-r}$ jest obrazem wektora $\boldsymbol{v}$ w odbiciu symetrycznym względem płaszczyzny stycznej do powierzchni i odpowiada pozbawionemu fizycznego sensu załamaniu światła, które przechodzi przez powierzchnię. Na rysunku F.3 widać (w górnej prawej części) niewidoczną z punktu $\boldsymbol{e}$ gałąź krzywej, której punkty, będące miejscami zerowymi funkcji $f$, odpowiadają takiemu "odblaskowi".
\\\\
\textbf{Uwaga:} Warstwice funkcji $f$ inne niż zerowa (tj. krzywe odpowiadające zbiorom rozwiązań równania $f(u,v)-e=d$ dla $c \neq 0)$ na ogół \textit{nie są} liniami odblasku (tj. obrazami prostych odbitych w powierzchni), choć mogą być stosowane do badania kształtu powierzchni. Ta sama uwaga dotyczy również opisanych dalej linii odblasku niezależnych od obserwatora.

W praktyce zwykle do oceny kształtu służy rodzina linii odblasku, odpowiadająca pewnemu zbiorowi prostych równoległych, leżących w jednej płaszczyźnie i równoodległych. Stosuje się je często w projektowaniu karoserii samochodu w ten sposób, że ustawia się go w pomieszczeniu, na którego ścianach są rozmieszczone świetlówki. Projektant (stylista) może ocenić estetykę kształtu linii odblasku. Zanim powstanie prototyp, który można umieścić w takim pomieszczeniu w celu ostatecznego sprawdzenia, projektant może poddać ocenie realistyczne obrazy karoserii z liniami odblasku, wykonane na podstawie komputerowej reprezentacji projektu, na przykład metodą śledzenia promieni.
\\***rys***\\
Na rysunku F.4 jest przedstawiona połowa powierzchni tylnej szyby samochodu osobowego, która składa się z 90 płatów wielomianowych stopnia $(3,3)$ (szare linie na rysunku przedstawiają brzegi tych płatów). Położenie obserwatora $\boldsymbol{e}$, które posłużyło do określenia narysowanych linii odblasku, pokrywa się ze środkiem rzutu perspektywicznego użytego do wykonania tego rysunku.

\textbf{Linie odblasku niezależne od obserwatora} są określone w sposób pokazany na rysunku F.5. Punkt $\boldsymbol{p}$ leży na takiej linii wtedy, gdy prosta prostopadła do powierzchni w tym punkcie przecina prostą $l$. Warunkiem należenia punktu $\boldsymbol{p}$ do linii odblasku jest liniowa zależność wektorów $\boldsymbol{n}$, $\boldsymbol{s}=\boldsymbol{a}-\boldsymbol{p}$ i $\boldsymbol{t}$ Zatem linia odblasku jest zbiorem miejsc zerowych funkcji 
\begin{equation}
f(u,v) = det[\boldsymbol{n}, \boldsymbol{s}, \boldsymbol{t}].
\end{equation}
Wektor $\boldsymbol{t}$ jest stały; funkcja opisująca wektor $\boldsymbol{s}$ dla płata klasy $G^n$ jest klasy $C^n$. Wektor normalny $\boldsymbol{n}$ i określony z jego pomocą wektor $\boldsymbol{r}$ są opisane przez odwzorowania klasy $C^{n-1}$, a zatem oba rodzaje linii odblasku na powierzchni klasy $G^n$ są klasy 
\\***rys***\\
$G^{n-1}$ (z wyjątkiem punktów, w których odpowiednia funkcja $f$ ma pochodne cząstkowe równe $0$). Podobnie jak izofoty, linie odblasku lepiej niż przekroje uwidoczniają nieciągłości kształtu. Na przykład, jeśli taka linia ma w pewnym pukcie załamanie (tj. nieciągły kierunek stycznej), to krzywizna powierzchni w tym punkcie jest nieciągła.

\subsection{Krzywizny powierzchni}
Funkcje opisujące krzywizny główne powierzchni, a także jej krzywiznę gaussowską i krzywiznę średnią, są stosowane w badaniu kształtu bardzo często. Obrazy tych funkcji są oczywiście niezależne od parametryzacji. Jeśli powierzchnia jest klasy $G^n$, to jej krzywizny są funkcjami $C^{n-2}$, a ich warstwice na powierzchni-krzywymi klasy $G^{n-2}$.
***rys***
(***)
Na rysunku \verb|\ref{fig:F.6}| jest pokazany obraz krzywizny gaussowskiej szyby (tej samej co na rysunku \verb|\ref{fig:F.4}|). Paleta użyta na obrazku z lewej strony jest ciągła i ułatwia ocenę rozkładu krzywizny na powierzchni. Obrazek z prawej strony, wykonany przy użyciu palety nieciągłej, między innymi lepiej uwidocznia obszar o ujemnej krzywiźnie gaussowskiej (obszar ten jest wypełniony kolorem czarnym i ciemnoszarym).
***rys***

Na rysunku \verb|\ref{fig:F.7}| w podobny sposób jest przedstawiona krzywizna średnia tej samej powierzchni. Można dowieść (zobacz np. [183]), że dla powierzchni klasy $G^1$ ciągłość krzywizny średniej jest warunkiem koniecznym i dostatecznym ciągłości $G^2$ tej powierzchni (ciągłość krzywizny gaussowskiej jest tylko warunkiem koniecznym ciągłości $G^2$). Dlatego obrazy krzywizny średniej, takie jak na rysunku, mają podstawowe znaczenie w projektowanie powierzchni o ciągłej krzywiźnie.

\section{Krzywe charakterystyczne}
Mając daną dowolną funkcję kształtu, możemy uzyskać obraz jej warstwic przez użycie palety nieciągłej. Metoda ta jest jednak cość czasochłonna i w wielu przypadkach nieodpowiednia. Na przykład każda z linii odblasku *** z rodziny określonej przez ustalony zbiór prostych jest warstwicą innej funkcji kształtu i trudno byłoby w ten sposób przedstawić na jednym rysunku wszystkie linie z takiej rodziny. Ponadto funkcje kształtu nie wyczerpują listy narzędzi wizualizacji powierzchni. Aby ocenić kształt powierzchni, określa się różne pola wektorowe, zależne od jej kształtu (czyli \textbf{wektorowe pola kształtu}). Obrazy krzywych otrzymanych w wyniku całkowania takich pól mogą być poddane ocenie, podobnie jak obrazy funkcji kształtu i ich warstwic.

\subsection{Całkowanie krzywych charakterystycznych}	
Rozważmy pole wektorowe $\boldsymbol{\upsilon}: A \to\mathbb{R}^2$ określone w pewnym obszarze $A$ zawartym w dziedzinie powierzchni \textbf{p}. Jeśli pole to jest różne od $\boldsymbol{0}$ w całym zbiorze $A$ i spełnia warunek Lipschitza, tj. jeśli istnieje stała $L$, taka że dla dowolnych dwóch punktów $\boldsymbol{b_1 i b_2} \in A$ zachodzi nierówność
\begin{equation*}\begin{tabular}{l}
$||\upsilon(b_1)-\upsilon(b_2)|| \le L||b_1-b_2||$,\\
$||\upsilon(b_1)-\upsilon(b_2)|| \le L||b_1-b_2||$,
\end{tabular}\end{equation*}
to dla każdego punktu tego zbioru istnieje przedział $(a, b)$ i przechodząca przez ten punkt regularna krzywa $\boldsymbol{b}: (a, b)\to A$, taka że dla każdego $t\in (a, b)$ zachodzi równość $\boldsymbol{b}'(t)= \boldsymbol{\upsilon(b}(t))$.

Znając krzywą $\boldsymbol{b}$, możemy narysować jej obraz $\boldsymbol{c}=\boldsymbol{p(b)}$ na powierzchni. Tak jak w przypadku funkcji kształtu i ich warstwic, dla oceny powierzchni nie ma znaczenia parametryzacja krzywej $\boldsymbol{c}$, tylko jej kształt. Dlatego podczas numerycznego wyznaczania krzywych całkowych pole $\boldsymbol{\upsilon}$ możemy dowolnie skalować, tj. mnożyć przez dowolne funkcje, pod warunkiem zachowania jego ciągłości. Często dokonuje się \textbf{unormowania} pola $\boldsymbol{\upsilon}$, czyli podzielenia każdego obliczonego wektora z tego pola przez jego długość. Krzywa całkowa pola unormowanego ma w dziedzinie płata prędkość jednostkową.

Niech $\upsilon(u,\upsilon)=[x(u,\upsilon), y(u,\upsilon)]^T$. Wtedy równanie linii pola $\boldsymbol{\upsilon}$ możemy przedstawić w postaci układu dwóch równań różniczkowych zwyczajnych:
\begin{equation}\label{eq:F.5}
\frac{d}{dt}u=x(u,\upsilon),\qquad\frac{d}{dt}\upsilon=y(u,\upsilon),
\end{equation}
Wybierając dowolny punkt $[u_0,\upsilon_0]^T\in A$ określamy dla układu \eqref{eq:F.5} zagadnienie początkowe $\boldsymbol{b}(0)=[u_0,\upsilon_0]^T$, które można rozwiązywać numerycznie. Najprostszy sposób rozwiązywania bierze się z rozwiązania przybliżonej równości $(\boldsymbol{b}(t+h)-\boldsymbol{b}(t))/h\approx\boldsymbol{b}'(t))=\boldsymbol{\upsilon(b}(t))$ ze względu na $\boldsymbol{b}(t + h)$ i jest znany jako \textbf{metoda Eulera}. Rozwiązując równania \eqref{eq:F.5} tą metodą, wybieramy \textbf{krok całkowania} $h$ i obliczamy współrzędne
\begin{equation*}
u_i=u_{i-1}+hx(u_{i-1},\upsilon_{i-1}),
\end{equation*}
\begin{equation*}
\upsilon_i=\upsilon_{i-1}+hy(u_{i-1},\upsilon_{i-1}),
\end{equation*}
kolejnych punktów $\boldsymbol{b}_i=[u_i,\upsilon_i]^T$, leżących na krzywej stanowiącej przybliżone rozwiązanie naszego układu równań różniczkowych. Rysunek krzywej otrzymamy, wyświetlając łamaną (dokładniej: rzut równoległy lub perspektywiczny łamanej, której wierzchołkami są punkty $\boldsymbol{p}(u_i,\upsilon_i))$.

Przyjmując $h>0$, a potem $h<0$, możemy otrzymać krzywą „po obu stronach” punktu $[u_i,\upsilon_i]^T$. Rozwiązanie jest tym dokładniejsze (mowa tu zarówno o dokładności metody całkowania, jak i o dokładności przybliżenia poszukiwanej krzywej przez łamaną), im mniejsza jest wartość bezwzględna parametru $h$ (ale ze zmniejszeniem $h$ rośnie czas obliczeń). Często metoda numerycznego całkowania zawiera procedurę wyboru długości kroku do wyznaczania kolejnych punktów. Całkowanie kończymy wtedy, gdy dojdziemy do brzegu dziedziny płata, trafimy na osobliwość (tj. punkt, w którym $\upsilon=0$), otrzymamy krzywą zamkniętą albo przetniemy inną krzywą całkową (co może być spowodowane błędami metody całkowania i błędami zaokrągleń). Niestety, z wyjątkiem pierwszej z wymienionych sytuacji, zaprogramowanie badania warunku stopu nie jest łatwe.

Regularna parametryzacja płata wielomianowego lub wymiernego jest jego parametryzacją o maksymalnej gładkości i na ogół nie sprawia kłopotów podczas całkowania krzywych. W przypadku powierzchni sklejanych (na przykład B-sklejanych albo złożonych z płatów połączonych z ciągłością $G^n$) ze szczególną uwagą należy zaprogramować „przekraczanie” wspólnego brzegu połączonych płatów. Właśnie w tych miejscach występują nieciągłości na przykład kierunku stycznej do krzywej (krzywa, która jest warstwicą funkcji kształtu może być w takim miejscu nawet nieciągła). Dlatego, dochodząc do takiego miejsca, należy dobrać ostatni krok całkowania tak, aby otrzymać punkt $\boldsymbol{b}_i$ na brzegu płatów (co sprowadza się do rozwiązania odpowiedniego równania nieliniowego). Otrzymany punkt na brzegu będzie punktem początkowym łuku krzywej leżącego na sąsiednim płacie. Takie postępowanie jest potrzebne do otrzymania obrazu, na którym nieciągłości kierunku stycznej nie ulegną zniekształceniu.

Zamiast metody Eulera w praktyce lepiej jest stosować inne metody. Wyróżniają się wśród nich \textbf{metody różnicowe} i \textbf{metody Rungego-Kutty}. Obie te klasy metod zawierają metodę Eulera jako przypadek szczególny. Odpowiednio dobrana metoda dla ustalonego kroku całkowania daje znacznie mniejsze błędy. Dokładne przedstawienie tych metod wykracza poza ramy tej książki, ale można je znaleźć w podręcznikach, np. [61], [166].

W wizualizacji kształtu powierzchni pewien kłopot sprawia wybór krzywych, które mają zostać wyznaczone i narysowane. Musi ich być dostatecznie dużo, aby projektant mógł ocenić każdy interesujący go fragment, jednak zbyt duże ich zagęszczenie zaciemnia obraz i utrudnia jego interpretację. Automatyczny wybór punktów, od których zaczyna się całkowanie, jest dość trudny, jeśli ma zapewnić w miarę równomierny rozkład krzywych na powierzchni. Inna możliwość to metoda „ręczna”; projektant może użyć np. myszy do wskazania interesujących go punktów.

\subsection{Warstwice i linie najszybszego spadku funkcji kształtu}
Niech $f:A\to \mathbb{R}$ oznacza dowolną funkcję (funkcję kształtu) klasy $C^1$ (dokładniej, funkcja powinna mieć ciągłe pochodne cząstkowe w dziedzinie płatów, np. wielomianowych, z których składa się dana powierzchnia). Aby wyznaczyć warstwicę tej funkcji przechodzącą przez punkt $[u_0, \upsilon_0]T$, można rozwiązać zagadnienie początkowe
\begin{equation}\label{eq:F.6}\begin{tabular}{ll}
$\frac{d}{dt}u=-f_\upsilon(u, \upsilon),$ & $\frac{d}{dt}\upsilon=f_u(u, \upsilon),$\\$u_0=u_0,$ & $\upsilon(0)=\upsilon_0.$
\end{tabular}\end{equation}
Istotnie, pole wektorowe, którego krzywa całkowa jest rozwiązaniem tego zagadnienia, jest w każdym punkcie prostopadłe do gradientu funkcji $f$.

Podczas numerycznego całkowania równań \eqref{eq:F.6} mamy możliwość poprawiania błędów. W tym celu wystarczy po wykonaniu każdego kroku (czyli po obliczeniu kolejnego punktu $[u_i, \upsilon _i]^T$ dokonać poprawki, rozwiązując równanie $f(u, \upsilon)-c=0$ (z parametrem $c=f(u_0, \upsilon_0)$). Możemy rozwiązać równanie (z niewiadomą $s$)
\begin{equation}\label{eq:F.7}
f(u_i+sf_u(u_i+\upsilon_i), \upsilon_i+sf_u(u_i, \upsilon_i))-c=0,
\end{equation}
na przykład którąś z metod opisanych w p.\verb|\eqref{E.4}|, a następnie podstawić

$\begin{bmatrix} u_i\\\upsilon_i\end{bmatrix}:=\begin{bmatrix} u_i+sf_u(u_i, \upsilon_i)\\\upsilon_i+sf_\upsilon(u_i, \upsilon _i)
\end{bmatrix}$.\\*
Jeśli błąd metody, spowodowany przyjęciem zbyt dużego kroku $h$, jest zbyt duży, to może się zdarzyć, że w przedziale np. $[-h, h]$ równanie \eqref{eq:F.7} nie ma rozwiązania. W takim przypadku możemy odrzucić ostatni punkt i przyjąć mniejszy krok całkowania. Jeśli to nie pomaga (krok trzeba zmniejszyć poniżej ustalonej granicznej wartości), to całkowanie należy zakończyć; całkowane pole wektorowe ma prawdopodobnie osobliwość, która powoduje nieciągłość warstwicy funkcji $f$ albo jej kierunku stycznej.

Oprócz warstwie do oceny powierzchni bywają używane również \textbf{linie najszybszego spadku} funkcji kształtu. Kierunek stycznej do takiej krzywej na powierzchni jest prostopadły do warstwicy funkcji w rozpatrywanym punkcie. Dla ustalonej funkcji kształtu $f$ pole wektorowe $\boldsymbol{\upsilon} = [-f_v,f_u]^T$ określa w każdym punkcie zbioru $A$ kierunek styczny do warstwicy funkcji $f$. Wektor pola $\boldsymbol{\omega}$, które wyznacza kierunki prostopadłe do warstwicy $f$ na powierzchni, jest prostopadły do $\boldsymbol{\upsilon}$ w sensie iloczynu skalarnego określonego przez pierwszą formę podstawową powierzchni p.\verb|\eqref{(C.2.2)|, Możemy go zatem określić wzorem

$\boldsymbol{\omega}=
\begin{bmatrix} g_{21}f_\upsilon-g_{22}f_u\\-g_{1}f_\upsilon+g_{12}f_u
\end{bmatrix}$.
\\*
\textbf{Ćwiczenie.} Sprawdź, że $\boldsymbol{I(\upsilon, \omega)} = 0$.

Zauważmy, że całkując pole w numerycznie, nie mamy możliwości korygowania błędów, tak jak podczas wyznaczania warstwie. Przykład obrazu linii prostopadłych do przekrojów powierzchni z rysunku $\verb|\eqref{F.1}|$ jest pokazany na rysunku \verb|\eqref{F.8}|.

\subsection{Linie krzywiznowe}

Równanie $(B-\chi G)\textbf{y} = \textbf{0}$, w którym litery $G$ i $B$ oznaczają macierze pierwszej i drugiej formy danego płata powierzchni,a $\chi$ - krzywiznę normalną tego płata w kierunku wektora $y$, ma dwa niezależne liniowo rozwiązania $\textbf{y} = \textbf{y}_1$ i $\textbf{y} = \textbf{y}_2$, które określają kierunki główne (zobacz p. C.2.4*). Odpowiadają im krzywizny główne $\chi_1$ i $\chi_2 $. Jeśli są one różne w każdym punkcie obszaru $A$ (zawartego w dziedzinie płata), to możemy określić dwa pola wektorowe w tym obszarze. Pierwsze z nich jet określone przez unormowane wektory \textbf{y} odpowiadające mniejszej, a drugie - większej krzywiźnie głównej. W punktach kulistych (takich ze $\chi_1 = \chi_2 $) kierunek każdego z pól jest nieokreślony.

\textbf{Linie krzywiznowe} powierzchni to krzywe otrzymane w wyniku całkowania pól wektorowych określonych wyżej. Ponieważ kierunki główne są wzajemnie prostopadłe, więc każdy punkt powierzchni, niebędący punktem kulistym, przechodzą dwie linie krzywiznowe, których styczne w tym punkcie są do siebie prostopadłe. Obrazy linii krzywiznowych są często stosowane do oceny kształtu powierzchni.

Linie krzywiznowe powierzchni o ciągłej krzywiźnie są klasy $G^1$, ale ciągłość kierunku stycznej wszystkich linii krzywiznowych na powierzchni nie jest (***) warunkiem dostatecznym ciągłości $G^2$ tej powierzchni (ćwiczenie: wskaż przykład powierzchni o nieciągłej krzywiźnie , której wszystkie linie krzywiznowe sa klasy $G^1$). Nieciągłości kierunku stycznej linii krzywiznowych (czyli nieciągłości kierunków głównych) powierzchni złożonej z płatów. W pobliżu punktu kulistego kierunek stycznej do linii krzywiznowej zwykle szybko się zmienia, co jest spowodowane brakiem ciągłości kierunków głównych w otoczeniu takiego punktu.